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^ i Collisionless plasma instabilities are fundamental in magnetic field genera- 

tion in astrophysical scenarios, but their role has been addressed in scenarios 
where velocity shear is absent. In this work we show that velocity shears must 
be considered when studying realistic astrophysical scenarios, since these trigger 
o ' the collisionless Kelvin-Helmholtz instability (KHI). We present the first self- 

consistent three-dimensional (3D) particle-in-cell (PIC) simulations of the KHI 
in conditions relevant for unmagnetized relativistic outfiows with velocity shear. 
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d \ such as active galactic nuclei (AGN) and gamma-ray bursts (GRBs). We show 

the generation of a strong large-scale DC magnetic field, which extends over the 
entire shear-surface, reaching thicknesses of a few tens of electron skin depths. 
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t^ I and persisting on time-scales much longer than the electron time scale. This DC 

magnetic field is not captured by MHD models since it arises from intrinsically 

^ \ kinetic effects. Our results indicate that the KHI can generate intense magnetic 

^ I fields yielding equipartition values up to e^/e^ c^ 10~^ — 10~^ in the electron 

time-scale. The KHI-induced magnetic fields have a characteristic structure that 
will lead to a distinct radiation signature, and can seed the turbulent dynamo 
amplification process. The dynamics of the KHI are relevant for non-thermal 

rS ■ radiation modeling and can also have a strong impact on the formation of rela- 

Cs3 I tivistic shocks in presence of velocity shears. 

Subject headings: magnetic field generation - plasma instabilities - Kelvin-Helmholtz 
instability 
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Introduction 



Magnetic fiel d generation in astr ophysical scenarios, such as AGN and GRBs, is not 
fully understood ( jColgate et al. II2OOII ). While such phenomena are of fundamental interest, 
they are also closely related to open questions su ch non-thermal radiation emission and cos- 
mic ray acceleration ( JBhattacharjee et al. II2OOOI ). Recently, collisionl ess plasma effects have 
been proposed as candidate me chanisms for magnetic field generation ( JGruzinov fc Waxman 



1999l : lMedvedev fc 



instability (IWeibel 



iOeb 



1959 



I999I). mediating; t he formation of relativist i c sho c ks via the Weibel 



Silvaet alJl2003h . The KHI f lD^Angelo I Il965l : iGruzinov I I2OO8I : 



Zhang et al. I l2009l ) should also be considered since it is capable of generating large-scale 
magnetic fields in the presence of strong velocity shears, which naturally originate in ener- 
getic matter outbursts in AGN and GRBs, and which are also present whenever conditions 
for the formation of relativistic shocks e xist. These larg:e-scale fields niay al so be further 



amplified by the magnetic-dynamo effect (JGruzinov 



2008: 



Zhang et al. Il2009h . 



Recent kinetic simulations have focused on magnetic field generation via electromag 
netic plasma instabilities in unmagnetized fiows without velocity shears. 3D PIC Simula 



tions of Weibel turbul ence ( ISilva et all l2003l : iFonseca et all l2003l : iFrederiksen et al. 1 12004 



Nishikawa et al. II2005I ) have demonstrated subequipartition magnetic field gen eration. The 
Weibel instability has been shown to be critical in mediating relativistic shocks ( jSpitkovsky 
20081 : Martins et al.ll2009l ). where a Fermi-like particle acceleration process has also been 
identified. These works have neglected the role of velocity shear in the fiow, which are 
an alt ernative raechan ism to generate sub-equipartition magnetic fields in relativistic out- 
fiows ( JGruzinov II2008I ). Furthermore, a shear fiow upstream of a shock can lead to density 
inhomoReneities via the KHI which may constitute important scattering sites for particle 
acceleration. In IZhang et al. I ( 120091 ). 3D magnetohydrodynamic (MHD) simulations of KH 
turbulence are discussed, observing magnetic field a mplification due to the KHI. In fact, the 
KHI setup is routinely used to test MHD models ( JMignone et al. I l2009l : iBeckwith et al. 
2OIII ). b ut the connection between the MHD description and the fully kinetic picture is still 
missing ( JGruzinov II2OIII ). However, the KHI contains intrinsically kinetic features, and so 
far 3D fully kinetic ab initio simulations have not been reported. 

In this Letter, we present the first self-consistent 3D PIC simulations of the KHI for both 
subrelativistic and relativistic scenarios of shearing unmagnetized electron-proton plasma 
clouds. We show that the KHI c ontains important kinetic features which are not captured 



in previous 3D MHD simula tions ( JKeppens et al. lll999l : IZhang et al. II2009I : iMignone et al. 



2OO9I : iBeckwith et al. II2OIII ). namely the transverse KHI dynamics, and the generation of a 
large-scale DC magnetic field at the shear region. This large-scale field can reach mC levels 
for typical parameters of the interaction of a relativistic fiow with the interstellar medium 



(ISM). Furthermore, our generalization of the KHI hnear theory ( JGruzinov II2008I ) to include 
arbitrary density jumps between the shearing flows allows us to conclude that the onset of the 
instability is robust to this asymmetry. The KHI can therefore operate in shears within the 
ejecta (similar density flows) and also between the eject a and the surrounding ISM (relative 
density ratios of 1 — 10), and it will likely operate at the same level (or stronger) than the 
Weibel instability, even for moderate velocity shears. In Section 2, we explore the density 
jumps effect on the behavior and features of the instability. The 3D PIC simulation results 
are presented and discussed in detail in Section 3. In Section 4, we discuss the saturation 
levels of the magnetic fleld, and conclusions are drawn in Section 5. 



2. Theoretical analysis 



The KHI linear theory, outlined in ( JGruzinov 1 120081 ) . is based on the relativistic fluid 
formalism of plasmas coupled with Maxwell's equations, and was analyzed for the particular 
case where the two shearing flows have equal densities. Realistic shears, however , are more 
likely to occur between diflFerent density flows ( JHardee et al. Ill992l : iKrause 1 120031 ). 



We have extended the analysis presented in ( JGruzinov 1 120081 ) for shearing electron- 
proton plasma flows with uniform densities n+ and n_, and counter- velocities +Vo and —Vq^ 
respectively. Here, the protons are considered free-streaming whereas the electron fluid 
quantities and flelds are linearly perturbed. The dispersion relation for electromagnetic 
waves is thus: 
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where u)' = uj/ujpj^ is the wave frequency normalized to the plasma frequency, with density 
n+, ujpj^ = (47rn+e^/7ome)-^/^. k^ = k c/cjp^ is the normalized wave number parallel to the 
flow direction, /3o = Vq/c and 70 = (1 — /^q)""^^^ is the Lorentz relativistic factor; e and rrie 
are the electron charge and mass, respectively, and c is the speed of light in vacuum. The 
unstable modes are surface waves, since their transverse wave number is evanescent. 

We now consider two scenarios. The similar-density (n+/n_ c^ 1) scenario, relevant for 
instance in internal shocks within the relativistic ejecta, where the densities of the shearing 
regions are comparable, and the n+/n_ > 1 scenario, relevant for external shocks associated 
with the interaction of the relativistic ejecta with the ISM. In the flrst (n+/n_ = 1), the 



modes with |/c'| < 1 are unstable, with a maximum growth rate F^^ 
for the fastest growing mode k[^^. 



Jcjp^ = i/Vs 



= ^max c/ujp^ = ^y 3/8/^0. The real part of uj^ vanishes 
over the range of unstable modes meaning that the unstable modes are purely growing 
waves. In the density-contrast (n+/n_ > 1) scenario, the shape of the spectrum is conserved 
for different density ratios, indicating that the general characteristics of the instability are 
unchanged. As n+/n_ increases, the bandwidth of unstable modes decreases, /c^^^ drifts to 
larger scales, and the growth rate T^^^^ decreases. Furthermore, the frequency uj^ acquires 
a real part over the range of the unstable modes, meaning that the growing perturbations 
propagate. The density-contrast unbalances the interaction between flows such that the 
growing perturbations are more strongly manifested in the less dense flow, and drift in the 
direction of this less dense flow. In the same-density case, the interaction between flows is 
balanced and the perturbations are purely growing. 

The normalized growth rate of the fastest growing mode, F'^^^^, decreases with the 
density contrast between flows. For n+/n_ ^ 1, the growth rate scales approximately as 
Tmax ^ (n+/n_)~-^/^ for both relativist ic and non-relativist ic shears, whereas in the high 
density contrasts (n+/n_ ^ 1), the growth rate scales as F(^^^ oc ( n^/n_)~-^/^ fo r non- 
relativistic shears (similar to the small cold beam-plasma instability, ( 10 'Neil Ill97ll )). and 
Tmax ^ (n+/n_)~-^/^ for highly relativist ic shears. These scalings show that the KHI will 
compete and operate on the same time-scales as the Weibel instability, and should thus be 
considered in realistic relativistic outflows where velocity shears are likely to be present. 



3D PIC simulations 



Numerical simulations were p erformed with QSIR 



and massively parallel PIC co de (IFonseca 



other astrophysical problems ( ISilva et al, 
2o3). 



et al 



200, 



i 



2003 



S, a f ully relativistic, electromagnetic, 
2008r>. which has been widely used in 



Medvedev et al. 



20051 . l2006l : 



Martins et al. 



We simulate 3D systems of shearing slabs of cold {vq ^ Vth^ where Vth is the thermal 
velocity) unmagnetized electron-proton plasmas with a realistic mass ratio rrip/me = 1836 
{nip is the proton mass), and evolve it until the electromagnetic energy saturates on the 
electron time-scale. In this Letter, we explore a subrelativistic (70 = 1.02) and a relativistic 
(70 = 3) KHI scenario, in the regime of two plasma slabs of equal density, wh ich is directl y 
relevant to internal shocks where 1 < 70 < 10 with comparable density shells ( jPiran II2005I ). 



The numerical simulations are prepared as follows. The shear flow initial condition is set 
by a velocity fleld with Vq pointing in the positive Xi direction, in the upper and lower quarters 



of the simulation box, and a symmetric velocity field with —Vq pointing in the negative Xi 
direction, in the middle-half of the box. Initially, the systems are charge and current neutral. 
In the subrelativistic case, the simulation box dimensions are 20 x 20 x 20 (c/cj^)^, where 
Up = (47rne^/me)'^/^ is the plasma frequency, and we use 20 cells per electron skin depth 
{c/(jjp). The simulation box dimensions for the relativistic scenario are 250 x 80 x 80 (c/ujp)^^ 
with a resolution of 4 cells per c/ujp. Periodic boundary conditions are imposed in every 
direction. 

The magnetic field structures generated in the subrelativistic and relativistic scenarios 
are displayed in Figures [1] and El respectively, and the magnetic field line topology is illus- 
trated in Figure [31 In the linear regime, the onset of the fiuid KHI, discussed in Section 
[21 occurs in the XiXs plane around the shear surfaces, generating a magnetic field normal 
to this plane, B2- This growing magnetic field component B2 is responsible for rolling-up 
the electrons to form the signature KH vortices (insets a2 and b2 of Figure [4]) . The typi- 
cal length-scale of the KHI modulations, observed in Figures [1] f and [2] f, agrees with the 
wavelength of the most unstable mode predicted by the fiuid theory (A^ax = 2 c/(jUp and 
Amax — 50 c/ujp in the subrelativistic and relativistic cases, respectively). The KHI modula- 
tions, however, are less noticeable in the relativistic regime because these are masked by a 
strong DC component with a magnitude higher than the AC component, which is negligible 
in the subrelativistic regime. This DC magnetic field {k' = 0), is not unstable according to 
the fiuid model as can be seen in equation ([T]). As the amplitude of the KHI modulations 
grow, the electrons from one fiow cross the shear surfaces and enter the counter-streaming 
fiow. Since the protons are unperturbed, due to their inertia, the current neutrality around 
the shear surfaces is unbalanced, forming DC current sheets, which point in the direction of 
the proton velocity. These DC current sheets induce a DC component in the magnetic field 
B2- The DC magnetic field is therefore dominant in the relativistic scenario since a higher 
DC current is setup by the crossing of electrons with a larger initial fiow velocity and also 
since the growth rate of the AC dynamics is lower by 7o compare d with the sub r elativ istic 



case. We stress that this DC field is not captured in MHD (e.g. IZhang et al. I ( 12009! )) or 
fiuid theories because it results from intrinsically kinetic phenomena. Furthermore, since 
the DC field is stronger than the AC field, a kinetic treatment is clearly required in order 
to fully capture the field structure generated in unmagnetized relativistic fiows with veloc- 
i ty shear. This chara cteristic field structure will also lead to a distinct radiation signature 
JMartins et al. Iboioh . 



Electr on density structures, which have not been reported in MHD simulations to our 



knowledge ( JKeppens et al. Ill999l : IZhang et al. II2009I : JMigno ne et al. II2009I : JBeckwith et al 



20111), emerge in the plane transverse to the fiow direction (insets al and bl of Figure [4]), and 



extend along the Xi direction forming electron current filaments. A harmonic perturbation in 



the B^ component of the magnetic field at the shear surfaces forces the electrons to bunch at 
the shear planes forming current filaments, which amplify the initial ma gnetic perturbation 



B:^. T his process is identical to the one underlying the Weibel instability ( JMedvedev fc Loeb 



19991 ) and leads to the formation of the observed transverse current filaments, along with the 
exponential amplification of B^ observed in Figure [H Figures [1] g and [2] g further show that 
the S3 magnetic field component shares a filamentary structure, underlining its connection 
in this process. The electrons undergoing this bunching process, slow down along their 
initial fiow direction. Again, since the protons are unperturbed at these time-scales, DC 
{kx^ = mode) current sheets are setup around the shear surfaces in a similar fashion to 
the longitudinal dynamics previously discussed. These current sheets induce a DC magnetic 
field in B2 (Figures [Hand [2] e), which is responsible for accelerating the evolving filaments 
across the shear surface, into the counter-propagating fiow. In the relativistic shear scenario, 
these filaments are strongly rotated due to the high intensity of S2, into the opposing fiow, 
leading to the formation of well defined finger-like density structures, as seen in inset bl 
of Figure [H These structures are less pronounced in the subrelativistic scenario due to the 
lower intensity of the DC component of B2 and to the slower Weibel-like electron bunching 
process (inset al of Figured]). Meanwhile, the current component J3, associated with the 
crossing motion of the electron current filaments along the X3 direction, induces the magnetic 
field component Bi (inset b of Figures [1] and [2]) . 

The growth rate of the magnetic field in the subrelativistic regime agrees with the the- 
oretical Fmax (Figured] a), whereas a significant mismatch is found in the relativistic regime 
(Figure d]b). This deviation is owed to the transverse dynamics of the KHI, which is not 
taken into account in the 2D theory. We performed 2D simulations matching the longitudinal 
{xiXs) and transverse (^2X3) planes of the 3D setups, in order to assess the independent evo- 
lution of the longitudinal and transverse dynamics of the KHI. The results of the 2D simula- 
tions of the longitudinal planes were in excellent agreement with the theoretical growth rates 
for both subrelativistic and relativistic cases, measuring 0.35 ujp and 0.07 cj^, respectively. 
In the 2D simulations of the transverse dynamics we measured T = 0.1 ujp and F = 0.3 ujp 
for the subrelativistic and relativistic scenarios, respectively. Therefore the full 3D evolu- 
tion of the KHI is mainly determined by the transverse dynamics in the relativistic regime, 
in contrast to the subrelativistic regime where the longitudinal dynamics is dominant. We 
stress that these growth rates are faster /co mparable to other collisionless plasma p rocesses 



that would occur in interpenetrating fiows ( JKazimura et al. Ill998l : ISilva et al.ll2003l ). 



4. Electron KHI saturation 

At later times, the growing KH perturbations begin to interact nonlinear ly, ultimately 
guiding the system into a turbulent state. Eventually, all large-scale shear surfaces in the 
electron structure become extinct, and the instability saturates. This stage is reached at 
roughly t c^ 100 /up in both subrelativistic and relativistic scenarios (Figured]). In this tur- 
bulent state, the drift velocity of the electrons vanishes in favor of heating and the magnetic 
field is mainly sustained by the proton current sheets close to the shear-surfaces. Here, most 
of the magnetic field energy is deposited in ^2, which has a uniform DC structure that 
extends throughout the entire shear-surfaces (which can be extremely large-scale in realistic 
shears), with a characteristic transverse thickness Lgat = (^ c/ujp^ where a typically ranges 
from 5 to 20, as measured in the simulations. The AC modulations in the magnetic field 
structure are at this time negligible compared to the DC component. Using Faraday's equa- 
tion and neglecting the displacement current term, we may estimate the maximum amplitude 
of the magnetic field as 5dc ^ SyrenoLsat'^o/c. In the case of relativistic shears, the estimate 
yields 

5dc - 1.6 a ^/no[cm-^] [mC], (2) 

where no is the plasma density. In the case of a relativistic plasma shear of ISM density, 
no — 1 cm~^, and assuming {a) ^ 10, the maximum amplitude of the DC magnetic field is 
on the order of 10 mC, spreading over the entire shear surface with a thickness of 5 x 10^ 
cm. Our simulations confirm this estimate, since they are performed in normalized units 
(time is normalized to cj"-^, space to c/ujp^ and magnetic field to rrieC ujp/e)^ showing similar 
Sdc levels. The 5dc estimate can also be used to determine the average equipartition value 
^b/^p (ratio of magnetic to initial particle kinetic energy) of the system as 

Cb 1 me 7o + 1 2 

Cp 8 rup % 

where tp = no{mp + me)c^{jo — 1)^x1^x2^x3 {Lx. being the simulation box dimension in 
the Xi direction), and r = Lsat/I/^s- We measure from the simulation at saturation a c^ 4 
(r c^ 0.2) for the subrelativistic case, yielding e^/e^ = 4 x 10~^, and a c^ 15 {r c^ 0.19) for 
the relativistic case, yielding e^/e^ = 1.2 x 10~^, which are comparable with the simulations 
results (Figure H]). We may also estimate the maximum value of equipartition which is 
found at the shear surfaces by averaging around the interaction region, which has a typical 
thickness on the order of Lgat : 

^b\ 1 me 7o + 1 2 



^pJmax 8^P 70 



yielding 2 x 10 ^ for the subrelativistic case and 7 x 10 ^ for the relativistic case. A higher 
efl&ciency of conversion of particle kinetic energy to magnetic fields is observed in the rel- 



ativistic case since the thickness of the proton currents sheets (Lgat) is much higher than 
in the subrelativistic case. Most of the energy in the system, however, is stiU contained in 
the ions which remain unperturbed at these time-scales. We expect the system to reach 
higher levels of equipartition once the protons undergo the proton-scale KHI, which would 
occur at roughly tproton-KHi — 100 {mijvne)^!^ juj^ o^ 4000/6^^. The long-lived large-scale 
DC magnetic field can thus be sustained up to the pro ton time-sca le which is roughly long 



enough for prompt GRB emission and early afterglow (jPiran Ill999l ) 



5. Conclusions 

In this Letter, we present the first self-consistent 3D PIC simulations of the KHI in 
unmagnetized electron-proton plasmas and analyze their evolution on the electron time- 
scale. Our results show that the multidimensional physics of the KHI is extremely rich and 
that kinetic effects play an important role, in particular, in the transverse dynamics of the 
KHI (which is dominant over the longitudinal KHI dynamics in relativistic shears), and in 
the generation of a strong large-scale DC magnetic field. 

The transverse dynamics of the KHI consists of a Weibel-like electron bunching process, 
leading to the formation of electron current filaments which are then accelerated across 
the shear surface, forming finger-like structures. At the electron saturation time-scale, the 
magnetic field has evolved to a large-scale DC field structure that extends over the entire 
shear-surface, reaching thicknesses of a few tens of electron skin depths, and persisting on 
time-scales much longer than the electron time-scale. This field structure is not captured 
by MHD or other fiuid models, and will lead to a distinct radiation signature. We measure 
maximum equipartition values of e^/e^ c^ 2 x 10~^ for the subrelativistic scenario, and 
cb/^p c^ 7 X 10~^ for the relativistic scenario. These equipartition values, which are typically 
treated as a fr ee parameter in radiation m odeling^:, match the values inferred from GRB 



afterglows by (jPanaitescu fc Kumar I l2002l ) . Moreover, the onset of the KHI is robust to 



density asymmetries making it ubiquitous in astrophysical settings. The KHI may operate 
in n+/n_ ^ 1 regimes, relevant in GRB internal shocks, and also in n^/ri- ^ 1 regimes, 
which are important in external shocks. Future work will address the impact of the KHI in 
the formation of relativistic shocks in the presence of velocity shears. 
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Fig. 1. — Magnetic field structures generated by shearing subrelativistic e~p^ flows with 
7o = 1.02 taken at time t = 4Q/ujp. The 3D visualizations (a), (d) and (g) correspond to the 
magnetic fleld components 5i, B2 and S3, respectively. The 2D slices of the magnetic fleld 
intensity (b), (e) and (h) are taken at the centre of the box Xi = 10 c/ujp^ and slices (c), (f) 
and (i) are taken at X2 = 10 c/ujp. 
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Fig. 2. — Magnetic field structures generated by shearing relativistic e~p'^ fiows with 70 = 3 
taken at time t = 69 /(jUp. The 3D visuahzations (a), (d) and (g) correspond to the magnetic 
field components 5i, B2 and ^3, respectively. The 2D slices of the magnetic field intensity 
(b), (e) and (h) are taken at the centre of the box Xi = 125 c/cjp, and slices (c), (f) and (i) 
are taken at X2 = 40 c/ujp. 
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Fig. 3. — Magnetic field lines generated in (a) the subrelativisitc scenario, and (b) the 
relativistic scenario, at time t = 100/up. 
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Fig. 4. — Evolution of the equipartition energy e^/e^ for a (a) subrelativistic and (b) rela- 
tivistic shear scenarios. The contribution of each magnetic field component is also depicted. 
The insets in each frame represent two-dimensional slices of the electron density at t = AQ/uOp 
and t = 69/(jUp for the respective case. The red (blue) color represents the electron density of 
the plasma that fiows in the positive (negative) Xi direction. Darker regions in the colormap 
indicate high electron density, whereas lighter regions indicate low electron density. Slices 
for insets (al), (a2), (bl) and (b2) were taken at the center of the simulation box; (al) and 
(bl) are transverse to the fiow direction, and slices (a2) and (b2) are longitudinal to the fiow 
direction. 



